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Abstract —We present a distributed resource allocation strategy 
to control an epidemic outbreak in a networked population based 
on a Distributed Alternating Direction Method of Multipliers 
(D-ADMM) algorithm. We consider a linearized Susceptible- 
Infected-Susceptible (SIS) epidemic spreading model in which 
agents in the network are able to allocate vaccination resources 
(for prevention) and antidotes (for treatment) in the presence 
of a contagion. We express our epidemic control condition as 
a spectral constraint involving the Perron-Frobenius eigenvalue, 
and formulate the resource allocation problem as a Geometric 
Program (GP). Next, we separate the network-wide optimization 
problem into subproblems optimally solved by each agent in a 
fully distributed way. We conclude the paper by illustrating per¬ 
formance of our solution framework with numerical simulations. 


I. Introduction 

In this paper, we consider the problem of controlling the 
outbreak of an epidemic process in a networked population 
in the absence of a social planner. We propose a distributed 
solution that enables local computation of optimal investment 
in vaccines and/or antidotes at each node, to respectively 
reduce its infection rate and increase its recovery rate; and 
contain the spread of an outbreak. Our distributed solution is 
based on an iterative algorithm - a Distributed Alternating 
Direction Method of Multipliers (D-ADMM) algorithm. The 
distributed solution framework presented here is based on a 
Geometric Programming formulation of the epidemic control 
problem, and its convex characterization first proposed in [1], 

The fairly recent outbreaks of the Middle Eastern Res¬ 
piratory Syndrome (MERS) virus [2] and Ebola virus [3] 
have rekindled interest in analysis and control of epidemic 
outbreaks in networked populations. The analysis of an SIS 
epidemic process in arbitrary (undirected) networks was first 
studied by Wang et al. [4] using a discrete-time model. Key in 
[4] was the establishment of an epidemic threshold based on 
the spectral radius of the network adjacency matrix. Ganesh 
et al. in [5] studied a continuous-time SIS epidemic model 
where the relation between the spectral radius of the network 
adjacency matrix and the speed of spreading was further 
corroborated. A similar result was proposed in [6]. 

In the event of an epidemic outbreak, of interest is how 
to optimally allocate vaccines and/or antidotes across the the 
population to control the spread of the epidemic. Sometimes, 
the focus is to hasten the rate at which the contagion is 
contained; or to minimize the associated cost of containing 
the outbreak. Significant research attention has been given to 
this problem in the last decade. In [7], the authors proposed 
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a greedy immunization strategy - immunizing acquaintances, 
and showed it outperforms a case where vaccinations are 
applied randomly across nodes in the network. In [8] the 
authors presented an analysis of immunization strategy using 
the PageRank vector of the network adjacency matrix. Wan 
et al. in [9] presented a method to design control strategies 
in undirected networks using eigenvalue sensitivity analysis. 
A convex framework to compute the optimal allocation of 
vaccines in undirected networks using Semidehnite program¬ 
ming (SDP) techniques was discussed in [10]. In [1], we 
generalized the convex formulation via Geometric Program¬ 
ming techniques to weighted, directed, strongly connected 
networks to compute the cost and speed optimal allocation 
of vaccines and/or antidotes in directed strongly (and not 
necessarily strongly) connected networks. 

In the absence of a social planner (or central coordinator), 
centralized computation of optimal investments in vaccines 
and/or antidotes as was done in [ ] poses a challenge. Fur¬ 
thermore, the network might be too large for a centralized 
optimization scheme. These motivate the need for an optimal, 
decentralized framework to achieve the network-wide objec¬ 
tive of controlling an outbreak. In this paper, we propose a 
fully distributed ADMM algorithm that enables local computa¬ 
tion of optimal investment in vaccines and/or antidotes at each 
node to control the spread of an outbreak in a directed, network 
comprising heterogeneous agents. While the discussion in this 
paper focuses on strongly connected networks, the results are 
easily amenable to networks that are not necessarily strongly 
connected. 

We formulate the epidemic control problem as a distributed 
resource allocation problem, in that agents in the network 
carry out local computations informed by interactions with 
neighboring agents to determine their optimal investment in 
vaccines and/or antidotes, and must cooperate to complete a 
global task. Our work is similar in spirit to [11] where tests 
for distributed control of positive systems were presented. A 
notable difference and major challenge between the particular 
problems addressed in [11] and this paper is that our problem 
in its natural form is a nonconvex optimization problem. A 
game-theoretic formulation has also been studied in [ 12] where 
different equilibria were analyzed for undirected networks, 
contrasting our distributed solution that applies to directed 
networks with non-identical agents. 

Distributed optimization methods exist in the literature and 
have applied to several classes of problems; see, for example 
[13, 14, 15], We employ a D-ADMM, a dual-based method, 
to solve our constrained resource allocation problem. Dual- 
based methods are usually used when the local optimization 
of each agent can be done efficiently. The D-ADMM, as we 
will illustrate, first decomposes our original resource allocation 
problem into two sub-problems. These sub-problems are then 
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sequentially solved in parallel by each agent, with the associ¬ 
ated dual variables updated at each iteration of the algorithm, 
allowing for a fully distributed implementation. In contrast 
to [16], where a D-ADMM algorithm was proposed for an 
unconstrained optimization problem on an undirected network, 
our problem in its original form is a constrained optimization 
problem on a directed network, where the need for satisfaction 
of a global constraint is necessary. 

The rest of the paper is organized as follows - In section 
II, we present the notation, spreading model considered, as 
well as state the resource allocation problem. Section III 
comprises a GP formulation of the problem, alongside its 
convex characterization. The proposed distributed solution 
is presented in Section IV; with results from our solution 
illustrated in Section V. A summary and concluding remarks 
follow in Section VI. 

II. Preliminaries, Model and Problem Statement 
A. Graph Theory 

We define a weighted graph as Q = (V,£, W), where V = 
{t>i,..., v n } is a set of n nodes, £ C V x V is a set of 
ordered pairs of nodes called edges, and the function W : 
£ —> R++ associates positive real weights to the edges in £. 
The node pair ( Vj,Vi ) form an undirected edge in Q . If Q is 
a directed network, the pair ( Vj,Vi ) is an oriented edge from 
node Vj to node v t . For an undirected graph Q, we define the 
neighborhood of node Vi as Ni = {j : (Vj,Vi) £ £}. When 
the graph Q, is directed, we define the in-neighborhood of 
node Vi (that is, the set of nodes with edges pointing towards 
Vi), as Nl n = {j : ( Vj,Vi) £ £}. A directed path from Vj, to 
t’ij in Q is an ordered set of vertices (u,;,. Vi 2 ,... such 

that \Vi s> Vi B+1 ) £ £ for s = 1,...,(. A directed graph Q is 
strongly connected if, for every pair of nodes u,, Vj £ V, there 
is a directed path from v, to Vj. 

We denote the adjacency matrix of a directed graph Q by 
A = [ <Hj ] = VV (v.j, v t ) if the edge ( Vj,Vi) £ £ points from 
Vj to Vi, and a,;j = 0 otherwise. Given an n x n matrix M, 
we denote by Vi (M ),..., v„ (M) and Ai (M ),..., A ra (M) 
the set of eigenvectors and corresponding eigenvalues of M, 
respectively, where we order them in decreasing order of 
their real parts, i.e., R(Ai) > K(A 2 ) > ... > R(A n ). We 
respectively call Ai (M) and vj (M) the dominant eigenvalue 
and eigenvector of M ; and denote by p(M), the spectral 
radius of M, which is the maximum modulus across all 
eigenvalue of M. Vectors are denoted using boldface letters 
and matrices using capital letters. The letter I denotes the 
identity matrix and R(,z) denotes the real part of z £ C. 

hi. Spreading Model 

The spreading model considered in this paper is a 
continuous-time heterogeneous (SIS) epidemic model and 
closely follows the development in [17]. In the heterogeneous 
model, we assume a networked Markov process in which 
nodes in the network can be in one of two states - susceptible 
or infected. Over time, based on each agent’s infection rate (3 ), 
recovery rate Si and interactions with neighboring agents and 
their states, the (infected or susceptible) state of each agent 


evolves. We assume the agents are non-identical in the sense 
that each agent has distinct infection rate (3 t and recovery rate 
Sp, and that these rates can be controlled by either injecting 
vaccines or antidotes. 

The evolution of the spreading model is as follows. At each 
time-step, the state of node vi is a binary random variable 
Xi(t) £ {0,1}. The states JQ(t) = 0 and X,(t) = 1 respec¬ 
tively indicate that node Vi is susceptible and infected. Let 
the vector of states be X(t) = (A'i(f),..., X n (t)) T . Suppose 
node Vi is in the susceptible state at time t. Its probability 
of switching to the infected state depends on its infection 
rate /3 j, the state of its neighbors {Xj(t), for j £ Af- n }. 
Mathematically, the probability of node Vi switching from 
susceptible to infected state can be expressed as 


Pr {Xi(t + At) = 1| X f (t) = 0, X(t)) = 

y aijpiXj (t) At + o(At), (1) 
je AT? n 

where At > 0 is considered an asymptotically small time 
interval. Similarly, suppose node v t is in the infected; the prob¬ 
ability of node Vi returning to the susceptible state depends on 
its recovery rate Si, as well as the states of its neighbors with 
incoming connections. More formally, this is given by 

Pr (Xi{t + At) = 0| Xi(t) = 1, X(t)) = Si At + o(At). (2) 


Given an n— agent network and two possible states each agent 
can be in, the Markov process has 2" states. The exponentially 
increasing state space of this model poses computational 
challenges and makes this model difficult to analyze for large 
networks. To overcome this challenge, we use a mean-field 
approximation of its dynamics [6], By applying mean-field 
theory to the Markov chain, the dynamics of infection spread 
can be approximated using a system of differential equations. 
Suppose we define Pi (t) = Pr (A, ; (t) = 1) = E (2Q (£)), i.e., 
the marginal probability of node V{ being infected at time 
t. Then, the Markov differential equation [6] for the state 
Xi ( t) = 1 is: 


dpi ( t) 
dt 


(1 - pi (£)) Pi ^2 aijPj (t) - SiPi (t). 
3 = i 


( 3 ) 


Since i = 1 ,,n, we can represent (3) more compactly by 
a system of nonlinear differential equation with dynamics 

^ = (BA - D)p(t) - P(t)BAp(t), ( 4 ) 

where p (t ) = (pi(t),... ,p n (t)) T , B = diag(/3 j), D = 
diag (Si), and P (t) = diag(pi (£)). This ODE presents an equi¬ 
librium point at p* = 0, called the disease-free equilibrium. 
Readers are referred to [6, 18] for a detailed treatment and 
presentation of the mean-field approximation noted above. 


Proposition 1 . Consider the nonlinear dynamical system in 
(4), where the network adjacency matrix A > 0, B = 
diag(pi) > 0 and D = diag(Si) > 0. Suppose the real part of 
the largest eigenvalue of BA — D satisfies 

R(Ai (BA-D)) < -e, 


( 5 ) 
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for e > 0, the disease-free equilibrium is globally exponen¬ 
tially stable; that is, ||p(t)||< ||p(0)|| A'e -et , for K > 0. 

Proof: See [10, 17] for proof; where we showed that 
the linear dynamical system p (t) = (BA — D)p(t) upper- 
bounds the dynamics in (3); thus, the spectral result in (5) is 
a sufficient condition for the mean-field approximation of (4) 
to be globally exponentially stable. ■ 

C. The Distributed Resource Allocation Problem 

We present the epidemic control problem as a resource 
allocation problem in which vaccines and antidotes are respec¬ 
tively applied at each node to reduce their infection rates and 
increase their recovery rates within feasible intervals 0 < /3^ < 
Pi < Pi and 0 < 6 t < Si < Si. The specific values of and 
Si depend on the amount of vaccines and antidotes allocated 
at node Vi . We assume that vaccines applied at node v, have 
an associated cost fi(Pi)', and antidotes have a cost gpSf). 
The function fi(Pi), assumed to be convex and monotonically 
decreasing with respect to pi V i = 1,..., n, represents the 
cost of modifying the infection rates of node Vi to some pi 
within the feasible interval to stabilize the spreading dynamics. 
Similarly, the function gpSf) is monotonically nondecreasing 
with respect to Si V i = 1,..., n, and represents the cost of 
modifying the recovery rates of node v, to some 6, within the 
feasible interval to contain the spread of an epidemic outbreak. 

We focus on a rate-constrained resource allocation problem, 
in which the objective is to find the cost-optimal, local alloca¬ 
tion of vaccines and antidotes to achieve a given network-wide 
exponential decay rate in the probability of infection. That 
is, given a desired infection decay rate e, the objective is to 
locally allocate resources to each node Vi such that ||p (t)\\ < 
||P (0)|| K exp (—ef), K > 0; where p = (pi,... ,p„) T , 
via local computations and interactions with neighbors. The 
problem is formally stated below: 

Problem 1. Given a directed network with associated adja¬ 
cency matrix A, node cost functions {fi(Pi), gi(Si)}™ =1 and 
bounds on the infection and curing rates 0 < P .< Pi < Pi, 
and 0 < Si < Si < Si respectively, and an exponential decay 
rate e > 0. Locally (via interaction between neighboring 
nodes), determine the cost-optimal distribution of vaccines and 
treatment resources to attain the desired decay rate. 


Mathematically, the objective is to: 

n 

minimize E fi(Pi) + gi(Si) (6) 

subject to R(Ai (BA — D)) < —e, (7) 

P, < Pi < Pi, ( 8 ) 

Si<Si<Si , i = l,...,n, (9) 


where (6) is the total investment across all agents, (7) con¬ 
strains the decay rate to e, and (8)-(9) maintain the infection 
and recovery rates in their feasible limits. Our goal is to solve 
(6)-(9) in a fully distributed manner, with the computation 
of investment per node fi(Pf) + gpSf), carried out locally. 
After presenting a convex characterization of the problem, we 
illustrate how it decomposes nicely for a distributed solution. 


III. The Resource Allocation Problem as a GP 

In this section, we formulate the resource allocation problem 
in (6)-(9) as a Geometric Program and present its convex char¬ 
acterization. Geometric Programs (GPs) are a class of nonlin¬ 
ear, nonconvex optimization problems that can be transformed 
into convex optimization problems and efficiently solved using 
interior-point methods, yielding globally optimal solutions 
[19], Building blocks of a GP are monomial and posynomial 
functions. Let the vector x = (x\,... ,x n ) £ ®++ denote 
n decision variables. A monomial functions (of the variables 
xi,... ,x n ) is a real-valued function /(x) = cxffx^ .. .x’fy , 
where c > 0 is the coefficient of the monomial and a,i £ M 
are exponents of the monomial. Posynomial functions are 
sums of monomials; that is, a function of the form /(x) = 
Yjk =l c k xl lk x% 2k ... x °^ k , where c k > 0, for k = 1,2,..., K 
and cijk £ R for j = 1 ,n and k = 1 ,,K. For example, 
while Xi + x\ and xix^t: are posynomials, xi — x\ is not 
a posynomial. A GP in standard form is one in which 
a posynomial function is minimized subject to posynomial 
upper bound inequality constraints and monomial equality 
constraints. More formally, it is of the form 

minimize /o (x) 

subject to /j(x) <1, i = 1,..., m, (10) 

fli(x) = 1, i = 

where f are posynomial functions and gi are monomials. 
GPs in standard form are not convex optimization problems, 
since posynomials are not convex functions. However, with a 
logarithmic change of variables and multiplicative constants: 
y, = log a;*, bi = log ci, b ik = log c ik and a logarithmic 
change of the functions’ values, we can transform (10) to the 
following equivalent problem in the variable y: 

K 0 

minimize ho(y) = log E exp( a o fc y + b ok ) 

k=1 

Ki 

subject to hi(y) = log E exp(af fc y + b ik ) <0, Vi 

k =1 

Qi(y) = aTy + bi = o, i = l,..., m. 

(ii) 

Problem (11) is convex and can be efficiently solved in 
polynomial time. (See [19] for a more detailed exposition on 
GPs). Hence, if our cost function is a posynomial (or more 
generally, convex in log scale), and if the constraints can be 
represented as monomial and posynomial functions, (6)-(9) 
can be formulated as a GP in convex form. 

Via the Perron-Frobenius lemma, from the theory of non¬ 
negative matrices, we express the spectral constraint (5) in 
an equivalent form to transform it to a set of posynomial 
constraints in the decision variables. 

Lemma 2. (Perron-Frobenius) Let M be a nonnegative, 
irreducible matrix. Then, the following statements about its 
spectral radius, p (M), hold: 

1) p (M) > 0 is a simple eigenvalue of M, 

2) Mu = p (M) u, for some u € M" + , and 

3) p(M) = inf {AgK:Mu< Au for u £ K" + }. 
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Recall that our model assumes the contact network is 
directed, strongly connected. Since the adjacency matrix asso¬ 
ciated with a strongly connected, directed graph is irreducible. 
Lemma 2 holds for the spectral radius of the adjacency matrix 
of any positively weighted, strongly connected digraph. A 
corollary of Lemma 2 is the following: 

Corollary 2. Let M be a nonnegative, irreducible matrix. 
Then, its eigenvalue with the largest real part, Ai (M), is real, 
simple, and equal to the spectral radius p (M) > 0. 

In [ I ], based on Propositions 1 and Lemma 2, the following re¬ 
sult established the formulation of (6)-(9) as a GP comprising 
monomial and posynomial functions. 

Theorem 3. ([!]) Given a strongly connected graph Q with 
adjacency matrix A = [Aij], posynomial cost functions 
{fi {Pi) 1 9i (^i)}j=i> bounds on the infection and recovery 
rates 0 < p. < Pi < /3j and 0 < <5, < Si < Si, 
i = l,... ,n, and a desired exponential decay rate £. Then, the 
optimal investment on vaccines and antidotes for node Vi to 
solve Problem (6)-(9) are fi (P *) and gt ^A + 1 — <5*^, where 

A = max {e, 5i for i = 1,..., n} and Pf,S* are the optimal 
solution for Pi and Si in the following GP: 


n 


minimize ^ f k {/3 k ) + g k ( t k ) 
fc=1 

(12) 

Pi Yli -1 A ij U j + $i U i 
subject to - . f - , -< 1, 

(13) 

(A + 1 - ej Ui 


(ti + S^j j ^A + 1^ <1, 

(14) 

A + l — <^<<5j<A + l — (5^, 

(15) 

< Pi < Pi, i = 1,- - .,n. 

(16) 

Proof: See Theorem 12 in [l] for proof. 

■ 

For simplicity of notation and ease of reading, we re-express 

(12)-(16) as 


minimize V' f t (Pi) + g.^Si) 

i=l 

(17) 

Pi Xu-l A ij U i + $i U i 

subject to - - - < 1, 

Ui 

(18) 

5 < Si < S 

(19) 

P<Pi<P, i= 1, ■ ■ ■, n; 

(20) 

where the auxiliary variables A. t., .S, introduced in 

Theorem 

3 to express the spectral constraint as a set of posynomial 
functions have been factored into the upper and lower bounds 


on Si in (19). The rest of the paper will focus on developing 
a distributed solution to (17)-(20). 


A. Separability of (17)-(20) 

To implement a distributed solution, first note that the cost 
function and constraint functions of the GP in (17)-(20) is 
separable per agent Vi. Each agent is able to locally solve the 


following problem: 


minimize ffpf) + g^Sf) 

(21) 

Pi Xu-i AijUj + SiUi 

subject to 3 <1, 

Ui 

(22) 

6< Si <6 

(23) 

P < Pi <P- 

(24) 


Though separable per agent, not all the decision variables are 
local. In particular, for agent Vi to minimize its cost function 
it needs the value Uj from nodes in its neighborhood set. 
The need for Uj in computing the optimum cost of node Vi 
is explicit in (22). To address this problem, we employ the 
Alternating Direction Method of Multipliers algorithm, which 
is well suited for such distributed optimization problems. 


IV. Distributed solution 

A. Alternating Direction Method of Multipliers (ADMM) 

The ADMM algorithm is a dual-based method for solving 
constrained optimization problems in which an augmented 
Lagrangian function is minimized with respect to the primal 
variables, and the dual variables are updated accordingly. 
Recent surveys on augmented Lagrangian methods and the 
ADMM algorithm can be found in monographs by Schizas 
et al., [20], Bertsekas [21] and a survey by Boyd et al. 
[22]; where illustrations and solutions to different optimization 
problems via the ADMM algorithm are presented. The ADMM 
algorithm works by decomposing the original optimization 
problem into subproblems that can be sequentially solved in 
parallel by each agent, allowing for distributed solutions to 
large-scale optimization problems. 

The standard ADMM solves the following problem 

minimize f{x) + g{z) 

x,z (25) 

subject to Ax + Dz = c, 

where the variables x £ R", z £ R m , matrices A £ 
R pxn , p g and c g The Augmented Lagrangian 

function for (25) is given by 

L p (x,z, A) = f(x)+g{z)-X T {Ax+Dz-c)+^\\Ax+Dz-c\\ 2 , 

(26) 

where A is the Lagrange multiplier associated with the con¬ 
straint Ax+Dz = c and p is a positive scalar. The update rules 
for the variables x, z and A in the ADMM implementation is 
given by 

x(k + 1) = arg min L p (x, z(k), A(k)) (27) 

X 

z(k + 1) = argminL p (a:(fc + 1 ),z, A(fc)) (28) 

A (k + 1) = A (k) — p (( Ax{k + 1) + Dz{k + 1) — c) (29) 

The updates in (27) - (29) is similar to those of dual descent 
algorithms (see [23] for instance), except that augmented 
Lagrangian is used and penalty parameter p is used as the 
step size in the dual updates. 
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B. Resource Allocation via D-ADMM 

Our goal is to present a distributed solution to the GP in 

(21) -(24). To make (21)-(24) amenable to a distributed solution 
via the ADMM, we introduce variables u, £ R n representing 
a local copy of the global variable u = (ui,...,u n ) T in 

(22) at each node ft. We also introduce an auxiliary variable 
z ij, that enables communication and enforces consensus in the 
values of u,; and Uj for all neighboring nodes (ft, v 3 ) £ £. We 
interpret the auxiliary variables z VJ as being associated with 
the edge (ft, ft) with the goal of enforcing consensus of the 
variables u, and u ; of its adjacent nodes ft and Vj. With these 
new variables, and consensus constraint, we can reformulate 
(21)-(24) as 


minimize /*(/%) + ft(ch) 


,. .. Pi £j=i A iM + s * u l . 
subject to - - --- < 1 , 


(30) 


II U * = X > 

j =i 


S< Si < s 

p<pi<p. 

Uj = z and = z ij, ( Vi,Vj ) £ £, 


where the scalar uj is the jth entry of the local estimate Uj at 
node Vi. The constraints Uj = and u 7 = z t] imply that for 
all pairs of agents (ft, ft ) that form an edge, the feasible set of 
(30) is such that Uj = u ,. If the network is strongly connected, 
these local consensus constraints imply that feasible solutions 
must satisfy Uj = u 3 for all, not necessarily neighboring, 
pairs of agents ft and v 3 . 1 Normalization of the vectors Uj 
for i = 1 ,..., n is to ensure that the local estimates u, have 
the same direction. 

To compute the augmented Lagrange of (30), let the dual 
variable oft and 7 ^ be associated with equality constraints 
Uj = z ij, and Uj = z* j, V ( Vi,Vj) £ £ respectively. In 
the iterations of the D-ADMM algorithm, the dual variable 
updates are: 

otij(k + 1) = otij(k) + |(u i(k) - u j(k)) V j £ N(i) (31) 

7 ij(k + 1) = 7 ij(k) + ^(uj(fc) - Uj(fc)) V j £ N(i). (32) 

Let Tj(/c + 1) = Pi(k + l),Si(k + l),Uj(fc + 1). Further, let 
M k ) — £jetv(i)( a u ( k ) + 7ft(*0) V Vi £ V. Given Uj(0) £ 
K" and </>j(0) = 0 , the iterative computations and updates are 
summarized in Algorithm 1. 


V. Numerical Simulations 

In this section, we illustrate performance of the D-ADMM 
Algorithm 1 on two strongly connected directed networks, 
and briefly discuss its convergence. As a proof of concept, 
and to show that functional correctness of Algorithm 1, we 

'in contrast to existing literature on Distributed ADMM algorithm [24] 
[25], where only consensus of local estimates is required, the distributed 
resource allocation problem we consider in addition to being constrained 
locally, requires consensus in local estimates. 


Algorithm 1 Distributed ADMM 

1: Given initial variables ft(0), ft(0) £ R, u,(0), <ft(0) £ R n for 
each agent Vi £ V. 

2: Set k = 1 

3: repeat 

4: For all v z £ V 

<j>i(k + 1) = <j>i(k) + p (ui(k) - Uj(k)) (33) 

jeN(i) 


r i(k + 1 ) = arg min fi{/3i) + gi{5i) + </>j Uj 

fii,8i ,uj 


subject to 


+ P H Ui_ 

J'6 N(i) 

Pi ELi A a u i + ^ u j 


U i(k) + u j(k) ,,2 
-o-II 2 


n w 


3 =1 

5<5i<5 

P<Pi<~P- 


5: Set k = k + 1 

6: until "_i ivjjJlu;(fc) — Uj(k)||< g, for g arbitrarily small. 


will show that the investment in vaccines and antidotes for all 
agents in the network converge to the solution obtained in the 
centralized case. Further, we will show convergence of the of 
the dual variables <f>j for all agents v,. 

The D-ADMM algorithm was used to solve Problem 1 on 
a 8 -node network with the following parameters: epidemic 
threshold r c = (1 — S)/p(A); ft = 2 t c , (3 = 0.2/3 and 
S = 0.9, S = 0.5. These parameters were chosen in such 
a way as to ensure the DFE is unstable in the absence of any 
investment in vaccines and/or antidotes, and stable otherwise. 
We normalize the investment in vaccines and antidotes using 
the following quasi-convex functions 

(34) 

The functions are such that for /3j = /?j, fi(P t ) = 0 and for 
Pi = P v fi{Pi) = 1. Similarly, for S. t =5i , ft; (5,) = 1 and 
for Si = Si, 9i{Sd = 0- The Lagrangian penalty parameter p 
was chosen to be 4. 

We illustrate the distributed solution on an 8 -node directed 
network, where the probability of a directed edge between two 

points is 0.32. The values of 5 = 0.025 and S = 0.750 and 

P = 0.1142 and p = 0.4393. Since a feasibility constraint is 
p(BA — D) < 1, the upper and lower bounds of Pi and 6, 
are such that 2 when /3j = /3 and Si = (3, across all agents in 
the network, p(BA — D) = 0.4600. Further, when Pi = p 
and Si = <5, p{BA — D) = 1.6. And when /3j = /3 and 

Si = S, p(BA — D) = 1.04. Finally, when /3j = P and 

Si = S, p(BA — D) = 1.02. These bounds ensure that all 

2 Specifically, the upper and lower bounds of ft and <5j were chosen as 
follows and based on the spectrum of A. We have <5 = 0.8, 5 = 3.9 X 3L, 
T c = (2/10)/Ai(A), 0 = 4t c , 0 = 30%ft all chosen in a way that 
ensures infeasibility of a solution when all agents are assigned the maximum 
or minimum possible infection or recovery rates. 


£T - Pi (i-?0 











6 


agents in the network are not easily allocated resources to 
yield the minimum possible infection rate or the maximum 
possible recovery rates in the network, since the epidemic 
control criterion will be violated. With the above parameters, 
the optimal solution to the centralized problem using the cost 
functions in (34) was 1.9731. As can be seen in Figure 1, we 
observe convergence of the total investment in the distributed 
solution to that obtained in the centralized solution. Figure 
2 illustrates the convergence of errors in the local estimates 
u, across all agents in the network. It shows that the local 
estimates u, of the global variable at node each node 
converges to those of their neighbors and the agents reach 
consensus in their estimates. As illustrated in Figure 3, the 
dual variables <pj for all nodes converge. 

We illustrate the solution on a fairly larger strongly con¬ 
nected network comprising 20 nodes. As was done earlier, the 
values 5 = 0.025, 8 = 0.750, /3 = 0.0641 and /3 = 0.2464 are 
chosen such that when fa = (3 and 8, = 8, across all agents in 
the network, p(BA — D ) = 0.3500. Further, when fa = f3 and 
Si = 8, p(BA — D) = 2.000. And when fa = /3 and Si = 8, 
p(BA — D) = 1.0750. Finally, when fa = /3 and Si = 8 , 
p(BA — D) = 1.2750. 

The total investment obtained from the centralized solution 
for the 20-node network was 4.8868. In Figure 4, we find 
that the total investment obtained from the distributed solution 
converges to that obtained in the centralized solution. We find 
that the local estimates u, of the nodes reach consensus as 
illustrated in Figure 5, where we show convergence of the 
errors in estimates to zero. Finally, the dual variables, do 
indeed, converge as illustrated in Figure 6 . 


A. On convergence of the D-ADMM Algorithm 1 

The ADMM algorithm applied to distributed optimiza¬ 
tion problems have typically considered distributed problems 
where the only constraint is consensus in the local estimates 
of the agents; for example [22, 24, 26], With just a consensus 
constraint and smoooth, differentiable convex cost, explicit 
computation of the local decision variables is possible. This 
allows for analytical expression of the optimal iterates, which 
enables convergence rate analysis of the algorithm. 

The problem considered in this paper, however, is a con¬ 
strained optimization problem where, in addition to the con¬ 
sensus constraint, each agent also has three local constraints 
to satisfy to achieve a feasible solution at each iteration of 
the algorithm as specified in (30). This informed the use of 
numerical solvers for computing the optimal local variables at 
each node, as done in Algorithm 1. 

Remark 1. It is known that ADMM algorithm converges when 
applied to convex problems [22, 24], The convex characteri¬ 
zation of our resource allocation problem via GP presented in 
Section III guarantees that the D-ADMM solution in Algorithm 
1 converges. The use of numerical solvers in computing the 
local optimal solution at each node (in line 4 of Algorithm 1) 
makes it difficult carrying out a convergence rate analysis. 


VI. Summary 

In this paper, we proposed a fully distributed solution 
to the problem of optimally allocating vaccine and antidote 
investment to control an epidemic outbreak in a networked 
population at a desired rate. The proposed solution was a 
D-ADMM algorithm, enabling each node to locally compute 
it’s optimum investment in vaccine and antidotes needed to 
globally contain the spread of an outbreak, via local exchange 
of information with its neighbors. In contrast to previous 
literature, our problem is a constrained optimization problem 
associated with a directed network comprising non-identical 
agents. Since numerical solvers are used to solve convex 
subproblems at each node, a convergence rate analysis of the 
D-ADMM algorithm is not presented. However, it is known 
that the ADMM algorithm converges for convex problems 
[22]; further, illustration of our results via simulations in 
Section V show that that the D-ADMM algorithm converges. 
The proposed distributed solution to the vaccine and antidote 
allocation problem for epidemic control presents a framework 
to contain outbreak in the absence of a central social planner. 
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